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Abstract 

The non-equilibrium dynamics of the first stage of the reheating process, that is 
dissipation via particle production is studied in scalar field theories in the unbroken 
and in the broken symmetry phase. We begin with a perturbative study to one loop 
and show explicitly that the mechanism of dissipation via particle production cannot 
be explained with a simple derivative term in the equation of motion. The dissipative 
contribution is non-local and there does not exist a local (Markovian) limit at zero 
temperature. Furthermore, we show that both an amplitude as well as a one-loop cal- 
culation present instabilities, requiring a non-perturbative resummation. Within the 
same approximations, we study an 0(2) linear sigma model that allows to study dissipa- 
tion by Goldstone bosons. We find infrared divergences that require non-perturbative 
resummation in order to understand the long-time dynamics. 

We obtain a perturbative Langevin equation that exhibits a generalized fluctuation- 
dissipation relation, with non-Markovian kernels and colored noise. 

We then study a Hartree approximation and clearly exhibit dissipative effects re- 
lated to the thresholds to particle production. The asymptotic dynamics depends on 
the coupling and initial conditions but does not seem to lead to exponential relaxation. 

The effect of dissipation by Goldstone bosons is studied non-perturbatively in the 
large N limit in an 0(N) theory. Dissipation produced by Goldstone bosons dramati- 
cally changes the picture of the phase transition. We find the remarkable result that 



1 



for "slow-roll" initial conditions (with the expectation value of the field initially near 
the origin) the final value of the expectation value of the scalar field is very close 
to its initial value. Most of the potential energy is transferred to the production of 
long-wavelength Goldstone bosons. We find that the minima of the effective action (as 
determined by the final value of the expectation value of the scalar field) depend on the 
initial conditions. This in fact points out that in the case of broken symmetry, in the 
large N limit there are many extrema of the effective action. In the Hartree and large 
N approximation dissipation occurs in a collisionless regime, a situation very similar 
to that of Landau damping. We provide extensive numerical analysis of the dynamics. 

1 Introduction 

It is well appreciated that the dynamics of dissipation in scalar field theories is of great 
importance in a variety of settings. One of the most interesting of these is that of reheating 
of the universe after an inflationary epoch has passed|]l], Recall that at the end of new 
or chaotic inflation[§], defined by when the slow-roll conditions for the so-called inflaton 
field[|, 0(t), fail to obtain, the inflaton begins to oscillate about its true ground state. 

Since the inflaton is coupled to lighter fields (fermions or other scalars) as the scalar field 
oscillates around the true minimum it decays into these other particles, and eventually these 
particles thermalize via collisions and relax to an equilibrium state at high temperature. It 
is usually stated that such coupling gives rise to a term of the form T<p{t), where F is the 
decay rate of the inflaton into the lighter fields|P, 0. This term, which is usually put in by 
hand as a phenomenological result of the existence of open decay channels, acts as a friction 
term in the equation of motion for the inflaton, and converts the inflaton energy density into 
that of the lighter particles during its oscillations. These decay products are then supposed 
to thermalize, completing the reheating process^, |^. 

This picture bears closer scrutiny. Two questions that must be addressed are: (i) how does 
a time reversible theory generate dissipative dynamics, and (ii) can the phenomenological 
term T(p{t) be derived from first principles? 

The answer to the first question is relatively well known. To generate dissipative dynamics 
from a time reversible theory requires that some "coarse-graining" of the field degrees of 
freedom be done. Essentially, one must trace out (in the functional integral) degrees of 
freedom other than the one that has been deemed important to the dynamics. The tracing 
operation turns a closed system (that of all the field modes) into an open one (that of the field 
modes of interest) where the traced out modes now become the "environment" which couples 
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to the remaining modes. Another important necessary ingredient are non-equihbrium initial 
conditions. 

As envisaged in the original scenarios, the process of reheating consists of two different 
stages. During the first one, potential and kinetic energy of the scalar field is dissipated by 
the process of particle production. The second stage involves the collisional thermalization 
of these particles reaching a final equilibrium state at some temperature. 

What we want to understand in this article is the dynamics of the first stage in an a 
priori fashion, i.e. how to incorporate the effect of quantum fluctuations into the dynamics 
of the time dependent order parameter, and the process of particle production. 

This is the first necessary step to fully understand the reheating process in its full com- 
plexity. 

What we do here is the following. Starting from a renormalizable self interacting scalar 
field theory we isolate the expectation value (which by translational invariance only depends 
on time) and generate its equation of motion taking quantum fluctuations into account. 
Since we are motivated by reheating after an inflationary epoch, we treat the dynamics at 
zero temperature since the temperature is supposed to have red-shifted all the way to zero 
during the inflationary stage. 

What makes our approach different from others that have been proposed (see below) is 
that we are not considering an "in-out" expectation value, such as would be generated by the 
usual method of constructing the effective action. Rather, we use a method that generates 
the equations for an "in-in" expectation value of the field operator. This ensures that our 
equations of motion for the field expectation value are both real and causal. This has been 
done by others, most notably by Calzetta and Hu|P] and PazQ. 

After obtaining the renormalized, causal and real equations of motion that determine the 
dynamics we provide an extensive numerical analysis of the equations. 

The results from our analysis are quite interesting and we summarize them here. The 
equations of motion are rather complicated integro-differential equations connecting the field 
expectation value to the fluctuations of the non-zero momentum modes. We begin by ob- 
taining the equations of motion in an amplitude expansion as well as in the loop expansion 
up to one loop. After recognizing the "dissipative terms" in these approximations, we find 
that these cannot be simply replaced by a term of the form T(p{t). 

In the case where there is no symmetry breaking, we generate the effective Langevin 
equation for the system, in the one-loop approximation, and we find a multiplicative, non- 
Markovian kernel for the dissipational term as well as a "colored" noise term. Thus this is 
very different than the simple friction term described above. In fact, it is easy to see that 
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there is no limit in which the dissipational dynamics we find can be described by a term 
proportional to the time derivative of the field. 

In the case in which a discrete symmetry is broken, so that there are no Goldstone modes, 
the one-loop equation can be linearized in the amplitude of oscillations about the non-zero 
expectation value of the field. Even in this linearized approximation, while a term involving 
the first time derivative of the field does appear, it is convolved with a non-trivial kernel, 
that again has no limit (at zero temperature) in which it becomes a local term. 

We show that the loop expansion or the amplitude expansion for the dynamics of the field 
expectation value must break down at long times. In the loop expansion, the order h term 
has an amplitude that grows in time due to resonance phenomena induced by the existence 
of a two particle threshold. Thus, it will eventually dominate over the tree amplitude, 
and perturbation theory will then break down. Dissipation can only be understood beyond 
perturbation theory. We then use the Hartree approximation in the single scalar field case 
and the large N approximation in an 0{N) symmetric theory to try to understand the effects 
of quantum fluctuations on the oscillations of the fleld expectation value as well as particle 
production due to these same oscillations. The Hartree approximation reveals dissipative 
effects due to particle production and open channels. However, we flnd that the damping 
is not exponential, and asymptotically, the fleld expectation value undergoes undamped 
oscillations with an amplitude that depends on the coupling and the initial conditions. The 
same features are found in the large N limit in the 0(N) model in the case of unbroken 
symmetry. Here we flnd that the asymptotic oscillatory behavior can be expressed in terms 
of elliptic functions. 

The dissipation in the system is caused by a Landau damping |jTO| type of process. This 
is a collisionless process which in our case is associated with the particle production arising 
from oscillations of the fleld expectation value via parametric ampliflcation. These quantum 
fluctuations then react back on the fleld expectation value, but not in phase with it, thus 
damping the oscillations. This behavior is similar to that found in the damping of strong 



electric flelds in a collisionless regime |Tl 



An unexpected and quite remarkable result comes from the analysis of the 0{N) model 
in large N when the symmetry is spontaneously broken. We flnd that the fleld expectation 
value oscillations are damped very quickly in this case, much more so than in the unbroken 
situation. This by itself is perhaps not so unusual, since the existence of Goldstone modes 
allows the fleld to lose energy by radiating a large number of soft Goldstone particles. What 
is unusual, however, is that there are some initial conditions, corresponding to "slow-roll" 
behavior for which the fleld expectation value relaxes (via dissipation) to a value that is 
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very close to the origin. This final value is a minimum of the effective action for which the 
expectation value of the scalar field has moved very shghtly from the initial value. Most of the 
potential energy is converted into long-wavelength Goldstone bosons. Thus the symmetry 
breaking phase transition is actually dramatically modified by the dissipational dynamics. 
These initial conditions are generically "slow-roll" conditions in which the scalar field begins 
very close to the top of the potential hill and the couplings are extremely weak (in our case 
about 10~^ or smaller). This is one of the most striking results of this article. 

We provide extensive numerical evidence for all our assertions; we can follow the evolution 
of the field expectation value, the quantum fluctuations and we also compute the number of 
particles created via the method of Bogoliubov coefficients. 

In the next section we provide an outline of the formalism we use to generate the equations 
for the expectation value of the scalar field that include the effects of the ffuctuations of 
the quantum ffuctuations. We then investigate the meaning of these equations within the 
one-loop and the amplitude expansion. In section III and IV we utilize non-perturbative 
approximations such as the Hartree approximation and the large approximation in the 
case of the 0{N) model to perform our numerical calculations. After this analysis, we 
reconcile the dissipative dynamics that we find with the concept of time reversal invariance. 

Section V contains our conclusions, in which we compare our results to other results on 
dissipation and reheating obtained in the literature. Two appendices are included which 
contain some of the technical details of our calculations. 
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2 Non-equilibrium field theory and equations of mo- 
tion 



The generalization of statistical mechanics techniques to the description of non-equilibrium 



processes in quantum field theory has been available for a long time[|l^, |13|, |I^ |T5|, [T6| but 
somehow has not yet been accepted as an integral part of the available tools to study field 
theory in extreme environments. We thus begin by presenting a somewhat pedagogical 
introduction to the subject for the non-practitioner. 

The non-equilibrium description of a system is determined by the time evolution of the 
density matrix that describes it. This time evolution (in the Schroedinger picture) is deter- 
mined by the quantum Liouville equation: 

zh^ = [H{t),p{t)] (1) 

A non-equilibrium situation arises whenever the Hamiltonian does not commute with the 
density matrix. Here we allow for an explicitly time dependent Hamiltonian, which might be 
the case if the system is in an external time dependent background, for example. 

The formal solution of the Liouville equation is: 

p{t) = U{t,to)p{to)U-\t,to) (2) 

with p{to) the density matrix at some initial time to that determines the initial condition for 
the evolution. Ensemble averages of operators and correlation functions are obtained as 

{0{t)) = Tr[oUit,to)p{to)U-\t,to)]/Trp{to) (3) 
{0{h)0{h)) = Tr[OU{h,t2)OU{t2,to)p{to)U-\h,to)]/Trp{to) (4) 

The expressions above have a very intuitive meaning. To compute the ensemble average, 
take the initial state (or density matrix), evolve it forward in time from to t, insert the 
operator in question and evolve the state backwards in time back to to- For the correlation 
function {ti > t2), evolve the initial state to t2, insert the operator, evolve it further to 
ti, insert the second operator, and finally evolve the state backwards to to. In most cases 
of interest the initial density matrix is either thermal or a pure state corresponding to the 
ground state of some initial Hamiltonian. In both cases the initial density matrix is of the 
form: 

p(to) = e-^^- (5) 
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the ground state of Hi can be projected out by taking the hmit f3 — > oo. The system will 
undergo non-equilibrium evolution whenever the time evolution operator does not commute 
with the initial Hamiltonian Hi. It is convenient to introduce a time dependent Hamiltonian 
H{t) such that H{t) = Hi for — oo < t < to and H(t) = Hevoiif) for t > to, where Hevoiif) is 
the evolution Hamiltonian that determines the dynamics. This corresponds to the assump- 
tion that the system has been in equilibrium up to and will evolve out of equilibrium 
thereafter. Real time equilibrium correlation functions can be obtained if H is constant in 
time for all times. 

Correlation functions and ensemble averages can be obtained by considering evolution 
forward in time (from to ~^ — oo) insertion of operators and backwards to the original 
time, and finally down the imaginary time axis to to — i(3 to account for the initial thermal 
condition. Taking the trace amounts to identifying the initial and final field configurations 
and performing a functional integral over this configuration. Thus one is led to consider path 
integrals in the complex time plane. 

Correlation functions or operator averages may be obtained as usual by coupling sources 
and taking functional derivatives with respect to them. 

Thus we are led to consider the following generating functional, in terms of time evolution 
operators in presence of sources 

Z[J+, J^] = Tr \U{T - ip, T- J'^)U{T, T'; J-)U{T\ T; J+)] (6) 

with T — oo; T' oo. The denominator in is simply Z[0, 0, 0] and may be obtained in 
a series expansion in the interaction by considering 2' [0,0, J^]. By inserting a complete set 
of field eigenstates between the time evolution operators, finally, the generating functional 
Z[J^, J~, J^] may be written as 

with the boundary conditions <I>+(T) = ^f^iT-if^i) = $ ; <l>+(r') = $-(T') = $2 ; $"(T) = 
$^(T) = $1. This may be recognized as a path integral along a contour in complex time. 
As usual the path integrals over the quadratic forms may be evaluated and one obtains the 
final result for the partition function 




exp \ i [ dt \Cinti-i5/5J^) - Cint{i5/5J' 




X 



exp <i I dtC 



int{-i5/5J^) \ exp jdti jdt2Jc{ti)Jc{t2)G,{ti,; 
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where Jc stands for the currents evaluated along the contour. The Gc's are the Green's 
functions on the contour |l7|, |l^, |^(the spatial arguments were suppressed). 

In the limit T — oo, the contributions from the terms in which one of the currents is 
or J~ and the other is vanish when computing correlation functions in which the external 
legs are at finite real time, as a consequence of the Riemann-Lebesgue lemma. For this real 
time correlation functions, there is no contribution from the terms that cancel between 
numerator and denominator. Finite temperature enters through the boundary conditions 
on the Green's functions. For the calculation of finite real time correlation functions, the 
generating functional simplifies to|^^ §] 

Z[J+, J-] = exp \i dt [Cint{-t6/6J^) - Cint{iS/6J-)] | x 

'-j^ dt^j^ dt2Ja{tl)Jb{t2)Gab{tlM)\ (9) 

with a, 6 = +, — . 

The Green's functions that enter in the integrals along the contours in equations (P, §) 
are given by (see above references) 

{T$(fi,ti)$(f2,t2)) = -iG++{nM;r2M) (10) 
(f<l>(fi,ti)$(f2,t2)) = -iG-{n,tr,r2M) (11) 

(<l>(fi,ti)<l>(f2,t2)) = 2G+-(fi,ti;f2,t2) (12) 

with T and T the time ordering and anti-time ordering symbols. The Green's functions above 
are written in terms of the homogeneous solutions to the quadratic form of the Lagrangian 
(free fields) G^ , G^ as 

G++(fi,ti;f2,t2) = G'>(fi,ti;r2,t2)e(ti-t2)+G'<(fi,ti;f2,t2)e(t2-ti) (13) 

G—{n,t^;f2,t2) = G'>(fi,ti;f2,t2)e(t2-ti)+G'<(fi,ti;f2,t2)e(ti-t2) (14) 

G+-(fi,ti;f2,t2) = -G<(fi,ti;f2,t2) (15) 

G-+(fi,ti;f2,t2) = -G>(fi,ti;f2,t2) = -G'<(f2,t2;ri,ti) (16) 

G>(ri,ti;r2,t2) = z($(fi,ti)$(f2,t2)) (17) 

G<{n,T;f2,t2) = G>{n,T-t(3,;f2,t2) (18) 

The condition (0) is recognized as the periodicity condition in imaginary time and is a 
result of considering an equilibrium situation for t < to- The functions G^ , G^, which are 
the homogeneous solutions of the quadratic form, with appropriate boundary conditions, 
will be constructed explicitly below. 
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The Feynman rules are the same as in regular perturbation theory, with the proviso 
that there are two interaction vertices corresponding to the ± branches, with opposite signs 
(arising from the complex conjugation of the unitary time evolution operator for evolution 
backwards in time). 

This formulation in terms of time evolution along a contour in complex time has been used 
many times in non-equilibrium statistical mechanics. There are many clear articles in the 



literature using this techniques to study real time correlation functions^, 17, 19, EG, 21, 22 



and effective actions out of equilibrium 0, |2^, |2^, ^ |. 

This formulation has already been used by some of us previously to study the dynamics 
of phase transitions p6[ . 



In our analysis we will take the /? — > cx) limit (zero temperature) following the argument 
provided in the introduction. We will also limit ourselves to a Minkowski space-time rather 
than expanding FRW space-time by considering phenomena whose time variation happens 
on scales much shorter than the expansion time of the universe, where H is the Hubble 
parameter. 

2.1 Equations of motion: 

The method briefly discussed above allows the derivation of the effective equations of motion 
for a coarse grained field. To focus the discussion, consider the situation of a scalar field 
theory with Lagrangian density 

£ = -d^^d''<^ - -m^^'^ - (19) 
2^24! 

The "coarse grained" field (order parameter) is defined as 

<^(t) = ($(x,t)) = ^^ (20) 

where we have used translational invariance. We write the field as $(x,t) = (f){t) + ilj{x,t) 
with ip{x,t) the fluctuations, obeying {ip^x.t)) = and consider (pit) as a background field 
in the evolution Hamiltonian. The non-equilibrium generating functional requires 




with £o[^^] the free field Lagrangian density of a field of mass m. We can now "integrate 
out" the fluctuations, thus obtaining the non-equilibrium effective action for the "coarse 
grained" background field. 
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The linear, cubic and quartic terms in ijj are treated as perturbations, while the terms 
(f)'^{t){ip'^)'^ may either be treated as perturbations, if one wants to generate a perturbative 
expansion in terms of the amplitude of the coarse grained field, or alternatively, they may 
be absorbed into a time dependent mass term for the fluctuations, if one wants to generate 
a loop expansion. We will now study both cases in detail. 



2.2 Amplitude expansion 
2.2.1 Discrete Symmetry 

We treat the term as perturbation, along with the linear, cubic and quartic terms. 

The conditions 

(V'^(x,t)) = (22) 
will give rise to the effective non-equilibrium equations of motion for the background field. 



This is the generalization of the "tadpole" method ||28|1 to non-equilibrium field theory. Since 
the mass of the fluctuations is the bare mass, the essential ingredient for the non-equilibrium 
propagators introduced above are the spatial Fourier transforms of the homogeneous solu- 
tions to the free field quadratic forms: 

^^(^1,^2) = -^e-^-'=(*i-*^) (23) 
G^ih,h) = -^e*-'=(*i-*2) (24) 
ujk = VP + m2. (25) 

To 0{(f)^) the diagrams are shown in figures (l.a-l.h). Since the (++) and the ( ) 

propagators are independent, the diagrams (l.a-l.h) finally lead to the equation of motion 



6"^ ' ' 2^^ ' J (27r)3 2ajfc 



Notice that the diagrams of figures (l.e-l.f) contain the one-loop two-particle threshold 
typical of the two particle to two particle scattering amplitude. The last term in ( p6|) is 
the real-time expression for this one-loop contribution. We can perform an integration by 
parts in the time integral, discarding the contribution at t = — oo by invoking an adiabatic 
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switching-on convergence factor to obtain 



f ^^'^ ^ M \ /■ ^^'^ cos[2uk{t - 1')] 



The fourth term in the above equation is recognized as a mass renormahzation and the 
fifth term as the couphng constant renormahzation. Using a Fourier expansion for (j){t) 
and its derivative, it is straightforward to see that after integration in the time variable, 
the remaining integration in k is ultraviolet finite. There are several noteworthy features 
of this expression that point to a more complicated description of dissipative processes in 
field theory. The first such feature is an expected one. The "dissipative" contribution, 
that is, the last term in the above equation, has a non-Markovian (i.e. memory-retaining) 
kernel. Secondly, the equation is non-linear in the amplitude of the coarse grained variable. 
These features are in striking contradiction with a simple r<p term in the equation of motion. 
Originally0 such a term was argued to arise if the scalar field in question is coupled (linearly) 
to other fields in the theory and truly speaking such a situation does not arise within the 
present context. However the issue of a memory-kernel (non-Markovian) is quite general[pU|, 



and a legitimate question to ask is: is there a Markovian (local) limit of this kernel? Such a 
limit would imply that 

)C{t - t') ^ DS{t - t') (28) 

with K,{t — t') being the non-local kernel present in the last term of (^). However, we find 
that at small {t — t') the kernel has typical logarithmic divergences. 

If a local approximation were valid, the "dissipative constant" D may be found by inte- 
grating the kernel in time. A straightforward calculation shows that infrared divergences give 
a divergent answer for such constant, clearly indicating that there is no Markovian limit for 



this kernel. In reference pj] a Markovian limit was argued to be available in the high temper- 
ature limit; since we are working at zero temperature the approximations invoked there do 
not apply. Thus, this lowest order calculation reveals two conclusive features that will per- 
sist in higher orders and even in non-perturbative calculations (see below): the "dissipative 
contribution " to the equations of motion obtained by integrating out the fluctuations are 
typically non-linear and, furthermore, they do not allow for a Markovian (local) description. 
We will postpone a numerical analysis of these equations until we study the full one-loop 
equations of motion below. 
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2.2.2 Broken symmetry 

In the case of unbroken symmetry, the above equation of motion does not admit a hneariza- 
tion of the dissipative contribution. However, in the case when the symmetry is sponta- 
neously broken such a linearization is possible 0. In this case let = — /i^ and let us look 
for small oscillations around one of the minima. To achieve this, write (pit) = ^6fi'^/X + 6{t). 
The mass term for the fluctuations now becomes and we will only keep the linear terms 
in 6 in (pT]). We now follow the steps leading to the equation of motion obtained before. 
That is, (^P^) = is imposed to one-loop order. After integration by parts in the time inte- 
gral (again with an adiabatic switching-on convergence factor) there appear several tadpole 
contributions. Those that are independent of 6 renormalize 0o corresponding to a shift of 
the position of the vacuum expectation value, while those that are linear in 6 renormalize 
the mass. We thus obtain in the linearized approximation: 



(2^ 2^ 



s+2,i6+^-^xni m I ^ r ^-:~r =^ (29) 



with Uk = a/A;^ + 2/1^. This expression is similar to that found by Paz|^. In that reference 
the short time behavior was analyzed. In order to solve this equation for all times, we must 
specify an initial condition. We will assume that for t < 6{t < 0) = 6i and 6{t < 0) = 
and that equation (P^ ) holds for t > OPD|. 

Under this assumption, the linearized equation can be solved via the Laplace transform. 
This solution corresponds to summing Dyson's series for the propagator with the one-loop 
contributions depicted in figure (1) (for the linearized case). The Laplace transform exhibits 
the two-particle cut[p, Q in the analytic continuation of the transform variable s. However 
the inverse Laplace transform proves to be extremely difficult to carry out. We will content 
ourselves at this stage with a perturbative solution to 0{\). Thus we write 

S{t) = Si cos[V2fxnt] + A5i(t) (30) 

The equation is then solved order by order in A. The corresponding equation for 6i{t) is 
now easily solved via the Laplace transform. The transform exhibits resonances at ±-\/2/i_R 
and the large time behavior is dominated by the secular term Si{t) «i tsm[\/2ixjit]. Thus 
the conclusion from the linearized approximation is that it is valid only for short times; 
the long time behavior will not be captured by this perturbative linearized approximation. 
This conclusion will be confirmed with a numerical analysis of the full one-loop case below. 



The analysis provided in references H, M, EBI in which a dissipative behavior was observed 



is not quite consistent, because their solution contains higher order terms in the coupling 
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that are not warranted in the approximation considered. Furthermore, even if the imphed 
resummation is accepted the solution found by these authors is only valid at very short times, 
corresponding to the region in the s-plane (Laplace transform variable) far away from the 
resonances. The physics of this effect is clear if analyzed in terms of the Fourier transform. At 
short times the contribution of large frequencies is substantial and the two-particle threshold 
is available for dissipation. However at long times only low frequencies have an important 
contribution, the two-particle threshold is not available and the resonances below threshold 
dominate giving a growing amplitude. 

Furthermore, even in this case in which there is an explicit linear velocity dependence in 
the "dissipative kernel" , the same argument presented before applies. There is no Markovian 
(local) limit in which this term becomes simply F^, despite the fact that S has a linear 
coupling to the fluctuations. 

2.2.3 Continuum symmetries and Goldstone bosons 

An important issue that we want to study in detail in this article is the process of dissipation 
by Goldstone bosons. As discussed in the introduction, dissipation via particle production 
becomes effective whenever the transferred energies become larger than multiparticle thresh- 
olds. To lowest order, as shown in the perturbative calculation of the previous section, the 
lowest threshold is the two-particle one, at an energy twice the mass of the particle. In 
the case of a spontaneously broken, continuous symmetry, there will be Goldstone bosons 
and the thresholds are at zero energy. In this case any amount of energy transfer may be 
dissipated by the Goldstone modes. 

This effect may be studied in the 0(2) linear sigma model with Lagrangian density 

We now write 

<y{x,t) = ^l^+6{t)+x{x,t) (32) 
and use the tadpole method to impose: 

(7r(x,i)) = ; {x{x,t)) = (33) 

Carrying out the same analysis in terms of Feynman diagrams to linear order in S{t) we find 
the following equation of motion (after integrating by parts invoking an adiabatic switching 
on convergence factor and absorbing the local terms in proper renormalizations) 

5 + 2i,l + 6i,IXr dt'5{t') \^K^{t - t') + hc^it - t')} = (34) 
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with ICy^it — t') the same kernel as in equation ( pQ] ) and 

^^(,_0^/|l«^.,„|M,*-01 (35) 

where M is an infrared cutoff introduced to define the integral. This infrared divergence is the 
result of the fact that the threshold for Goldstone bosons is at zero energy-momentum. This 
result clearly reflects many important features of "dissipation" via Goldstone bosons. First, 
as before, the dissipative kernel cannot be described in a local (Markovian) approximation 
even in this linearized theory in which the arguments leading to a T6 term would be valid. 
Second and perhaps more important, the long time behavior is clearly beyond perturbation 
theory as the contribution from the non-local kernels is unbounded as a result of infrared 
divergences associated with Goldstone bosons. 



2.3 One-loop equations 



The full one-loop equations of motion are obtained by absorbing the terms (f)'^{t){ip^y in 



eg . (|2T|) as a time dependent mass term for the fluctuating fields. In order to determine the 
dynamics for the background field, we will assume that 0(t < 0) = 0o and that at t = 
the background field is "released" with zero velocity. Now we must find the corresponding 
non-equilibrium Green's functions. Consider the following homogeneous solutions of the 
quadratic form for the fluctuations 







(36) 
(37) 



with 



(38) 



The boundary conditions (|37D correspond to positive and negative U~ frequency modes 
for t < (the Wronskian of these solutions is 2iul). Notice that U^(t) = [U^(t)]*. In terms 
of these mode functions we obtain 



G^{ti, t2 



2u 



-,Ut{t,)U,{t2) 



(39) 



Now the equation of motion to one-loop is obtained from the diagrams in figures (l.a), (l.d) 
(but with the modified propagator including the contribution of 0(t)). Thus to one- loop we 
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find the following equations 



f^.+(0) = 1 ; 








(40) 

(41) 
(42) 



where we have restored the h to make the quantum corrections explicit. This set of equations 
clearly shows how the expectation value (coarse grained variable) "transfers energy" to the 
mode functions via a time dependent frequency, which then in turn modify the equations of 
motion for the expectation value. The equation for the mode functions, (^T|) may be solved 
in a perturbative expansion in terms of A(^^(t) involving the retarded Green's function. To 
first order in that expansion one recovers eq. (|26|) . 

Before attempting a numerical solution of the above equations it is important to un- 
derstand the renormalization aspects. For this we need the large k behaviour of the mode 
functions which is obtained via a WKB expansion as in references []26| , p7| and to which the 
reader is referred to for details. We obtain 

d'k |f/+(t)p 



(27r) 



A2 
8^ 



1 
8^ 



In 



+ finite 



(43) 



where A is an ultraviolet cutoff and k an arbitrary renormalization scale. From the above 
expression it is clear how the mass and coupling constant are renormalized. It proves more 
convenient to subtract the one-loop contribution at t = and absorb a finite renormalization 
in the mass, finally obtaining the renormalized equation of motion 



0(t)+m^0(t) + ^03(t) + ^0(t) 



6 



\unt)\' 



xlh 



327r2 



0(t)(02(t)-02(O))ln[A/«:]=O 

(44) 

In the equations for the mode functions the mass and coupling may be replaced by the 
renormalized quantities to this order. One would be tempted to pursue a numerical solution 
of these coupled equations. However doing so would not be consistent, since these equations 
were obtained only to order h and a naive numerical solution will produce higher powers of 
h that are not justified. 

Within the spirit of the loop expansion we must be consistent and only keep terms of 
order Ti. First we introduce dimensionless variables 



r/(t) 



A 



R 



rriRt ; q 



k 

rriR 



9 



87r2 



(45) 
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and expand the field in terms of g as 



v{t) = Vd{r) + gr]i{T) + ■ ■ ■ 
Now tlie equations of motion consistent up to 0{h) become 



Vdir) + Vdir) + riiiir) 







^i(^) + ^i(^) + 3^/cz(^)m(7-) + Vdir) 



A/mn 



q^dq 



[q^ + l + 3v!M 



1/2 



+ 



(46) 



(47) 



(48) 



Tfie solution to equation (^) is an elliptic function. The equations for the mode functions 
become 

+q' + l + SvUT) 



Kir) 







(49) 



with the boundary conditions as in eq.(^2[) in terms of the dimensionless frequencies and 
where for simplicity, we have chosen the renormalization scale k = ttir. The chosen boundary 
conditions r]{Qi) = rjQ ; r]{0) = can be implemented as 



Vd{0)=r]o ; rjci{0) = ; = ; rjiiO) = 



(50) 



In fig. (2) we show ?7i(t) with the above boundary conditions with r/Q = 1 and g = 0.1. A 
cutoff A/m/j = 100 was chosen but no cut-off sensitivity was detected by varying the cutoff 
by a factor 3. Notice that the amplitude grows as a function of time. This phenomenon can 
be understood as follows. A first hint was obtained in the case of the linearized equations 
for a broken discrete symmetry. There we learned that because of resonances below the two- 
particle threshold, the amplitude grows at long times. This is the behavior shown in fig. (2). 
An alternative and perhaps more convincing argument is the following. The mode functions 
(r) obey a Schroedinger-like equation with a potential that is a periodic function of time 
because the classical solution is periodic. Let us call the period of the classical solution T. 
Floquet's theorem]^ guarantees the existence of solutions that obey 



Kir + T) 



(51) 



The Floquet indices fi are functions of the parameters of the potential. The classical solution 
is an elliptic function, and in this case the Schroedinger equation for the modes may be shown 
to be a Lame equation with n = 2 (see ref.[Q) whose solutions are Weierstrass functions. 
That is, a two- zone potential, with two forbidden and two allowed bands for g > 0. The 
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Floquet indices /i are pure imaginary for large g, but they are real and positive for q near 
zero. That is, the long wavelengths belong to an unstable band. 

A more intuitive understanding at a simpler level ensues if we look at small oscillations 
near the origin for the classical solution rjciir) ^ ?7d(0) cos(r). Then the Schroedinger equa- 
tion for the modes becomes a Mathieu equation^l|], for which the dependence of the Floquet 
indices on the parameters {q ; i]ci{0)) is known ||3ll|. There are unstable bands in which the 
Floquet index has positive real part for certain values of these parameters. 

Since the one-loop correction involves an integral over all wave-vectors, the values of q in 
these bands give growing contributions to the one-loop integral. 

These instabilities of the one-loop equations preclude a perturbative analysis of the pro- 
cess of dissipation. 



2.3.1 Goldstone bosons at one-loop 

It proves illuminating to study the one-loop contribution to the equations of motion for the 
scalar field expectation value from Goldstone bosons. We proceed as in the previous case 
but now with the Lagrangian density of the 0(2) linear sigma model (|3TD. We now write 



cr(x, t) 
(7r(x,t)) 



; {xix,t))=0 



(52) 
(53) 



where the terms x^ctq ; Tr^cTg are now included with the mass terms for the respective fields. 
The procedure is exactly the same as in the previous case, but we now use fiR as the mass 
scale and introduce the dimensionless variables of eq. (|45|) in terms of this scale. Performing 
the proper renormalizations and a subtraction at t = we finally find the following equations 
for this case 



dr 



2V -V + V + 9V 



1 



+ -ln[A//iH](^'(r)-r^^(0)) =0 (54) 



The mode functions satisfy the equations 

+g2-l + 3r/2(r) 








(55) 
(56) 
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where the boundary conditions are the same as in ( ^Tj ) but in terms of the initial frequencies 



W^{q) ; WT,{q) . As argued previously, for consistency we have to expand the field as in eq. (^) 
and keeping only rjci in the equation for the modes. Notice that whereas the minimum of the 
tree level potential corresponds to 77^; = 1 and the stable region for the U modes is rj^i > 1/3, 
the stable region for the V modes (pion field) is for rj'^i > 1. These instabilities require the 
choice of special initial conditions, in particular for the initial frequencies. We choose 



W^iq) = Jq'^ + l + 3r]^iO) (57) 



W^{q) = ^q^ + l+7]^{0) (58) 

This choice corresponds to a gaussian initial state centered at ri{0) but with positive fre- 
quencies at 77(0) = 0. This state evolves in a broken symmetry potential from the position 
determined by ?7(0). Alternatively, one can think of this situation as changing the sign of 
the mass at time t = from positive to negative. For small oscillations of rjd around the 
minimum at rjmin = 1, every time that the classical field oscillates to the left of the minimum 
the "pion" field becomes unstable and grows, and thus its contribution to the fluctuations 
becomes large. This behavior is depicted in fig. (3) which shows ?7i(t) with the boundary 
conditions as in (^) and the initial frequencies (|57| , [58| ) for the mode functions where chosen 



with 77d(0) = 0-6 and A/fiji = 100 but no cutoff sensitivity was detected by increasing the 
cutoff by a factor of 2. This situation clearly exhibits Floquet indices with a positive real 
part, and i]ci{t) is periodic for small oscillations around rjd = 1. The unstable wave vectors 
form a band when rjci oscillates around this value. In each period of rjd the amplitude grows 
because of these unstable modes. 

This figure clearly shows a dramatic growth in the amplitude of the quantum correction 
at long times as a consequence of the instabilities associated with the Goldstone mode. As 
a consequence of this, perturbation theory must fail at long times. 

2.3.2 The Langevin equation 

A fundamental aspect of dissipation is that of decoherence which plays an important role 
in studies of quantum cosmology ||25||. Furthermore dissipation and fluctuations are re- 



lated by the fluctuation-dissipation theorem. This point has been stressed by Hu and 
collaborators The relation between dissipative kernels and fluctuations and correlations 
of bath degrees of freedom is best captured in a Langevin equation. 

In this section we present a derivation of the corresponding Langevin equation to one-loop 
order in an amplitude expansion. The flrst step towards deriving the Langevin equation is to 
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determine the "system" and "bath" variables, once this is done, one integrates out the "bath" 
variables obtaining a (non-local) influence functional p5|. 



3J] for the system variables. 
Starting with the generating functional of eq.(|^) (for the case of zero temperature (3 oo) 
we separate the zero mode (0^) from the field as 



d^xi)^{x,t) = 



(59) 



and consider (p as the "system" and as the "bath" to be integrated out. Because of the 
condition in eg. (^9]) , there are no terms linear in ijj^ in the expansion of the action in terms 
of these fields. Thus 



+ 0((^^)^(^ 



(60) 



where Q is the spatial volume, S the action per unit spatial volume, 5*0 the free field action 
and the terms (9(('?/'^)^, (ip^)^) will contribute at the two-loop level and beyond. To this 
order, the coupling between "bath" and "system" is similar to the bi-quadratic coupling 
considered in reference (see also[25|). Integrating out the ip^ fields in a consistent loop 
expansion gives rise to the influence functional |25|, |33| , Q for the zero modes. The one loop 



diagrams contributing to this functional up to (9((0^)^) are shown in fig. (4). In order to 
obtain the Langevin equation it is convenient to introduce the center of mass {(p{t)) and 
relative {R{t)) coordinates (these are the coordinates used in the Wigner transform of the 



coordinate density matrix) asp2|, p3, p5 



0±(t)=0(t)± 



R{t) 



Using the Green's functions in eqs. (|TO|-|T^) with (^ 
we find the effective action per unit spatial volume 



(61) 

and being patient with the algebra 
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s, 



effi 



R] = r dt\c[(l) + R/2]-C[ 

J — CO I 



+ -R{tm) 



dt'(f)\t" 



+ i- 



dt 



oo 
oo 



- R/2] - -Rirnt) 
d^k sm[2cjfc(t - t')] 



d^k 1 



dt'R{t)R{t')^{t)(f){t') 



oo 



d^k cos[2cUfc(t-t')] 
(2^ 2^1 



+ 0{R^]R^ 



(62) 



The higher order terms 0{R^\ ...) receive contributions from two and higher loops and give 
higher order corrections to the lowest order (one-loop) Langevin equation. The imaginary 
part of the effective action above (last non-local term) gives a contribution to the path 
integral that may be written in terms of a stochastic field as 
1 



exp 



dt / dt'R{t)K{t,t')R{t') 



oc / V^V[i\ exp 



dt^{t)R{t) 



V[^] = exp 



with 



-If 

lC{t,t') 



/oo 
dt'^{t)}C-\t,t')^it') 
'OO 

d^k cos[2uJkit - 1')] 



4 



■mm 



{2n) 



The non-equilibrium path integral now becomes (keeping track of volume factors) 



Z (X J V(l)VRViV[i\ exp jifi 



5, 



■reffi 



R] + / dt^{t)R{t) 



(63) 
(64) 

(65) 
(66) 



with Sreff[4>,R\ the real part of the effective action in ( p^ ) and with V[^] the gaussian 
probability distribution for the stochastic noise variable given in (|63D. 



The Langevin equation is obtained via the saddle point condition ||2^, ^ p5 

6Sreff 



6R{t) 



\R=0 



(67) 



leading to 



d^k 1 



(pit) I dt'(f)\t'^ 



d^k sm[2ujk{t - t') 
(2^~ 



m (68) 



where the stochastic noise variable ^(t) has gaussian correlations 



« m »= ; « mm »= nm = -^mm 



d^k cos[2LJk{t - t')] 



(69) 
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Here, the double brackets stand for averages with respect to the gaussian probability dis- 
tribution V[£]. We can see that the noise is colored (not delta function correlated) and 
multiplicative. By integrating by parts the "dissipative kernel" in eq. ( |BHD (the last non-local 
term) in the same way as done in eg. (pTj) , we can clearly see that the resulting "dissipa- 
tive kernel" and the noise correlation function obey a generalized fluctuation-dissipation 
theorem p5t In the broken symmetry case, in the linearized approximation (around the 
tree level minimum) and if the k-integral could be replaced by a delta function, we would 
obtain the usual fluctuation-dissipation relation. By taking the average over the stochastic 
noise of (|^) with the gaussian probability distribution (|^) one obtains the equation of mo- 
tion ( p7D by replacing the average of product of fields by the product of the averages (thus 
considering the field classical background). 

The higher order terms in the effective action (influence functional) give rise to modifi- 
cations to the noise correlations, making them non-gaussian and involving more powers of 
(f) in the kernels, in principle these corrections may be computed systematically in a loop 
expansion. 

Although this Langevin equation clearly exhibits the generalized fluctuation-dissipation 
theorem connecting the "dissipative" kernel to the correlations of the stochastic noise, it 
is a hopeless tool for any evaluation of the dynamics. The long range kernels and the 
multiplicative nature of the noise prevent this Langevin equation from becoming a useful 
tool. Its importance resides at the fundamental level in that it provides a direct link between 
fluctuation and dissipation including all the memory effects and multiplicative aspects of 
the noise correlation functions. This last correlation function is related to the decoherence 
functional p5[|. 



2.3.3 Failure of perturbation theory to describe dissipation 

This section has been devoted to a perturbative analysis of the "dissipative aspects" of the 
equation of motion for the scalar field. Perturbation theory has been carried out as an 
amplitude expansion and also up to 0{h), both in the broken and the unbroken symmetry 
case. In both cases we found both analytically and numerically that the amplitude of the 
quantum corrections grow as a function of time and that the long time behavior cannot be 
captured in perturbation theory. This failure of perturbation theory to describe dissipation 
is clearly understood from a very elementary but yet illuminating example: the damped 
harmonic oscillator. Consider a damped harmonic oscillator 

q + Tq + q = (70) 
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with r ~ 0{\) where A is a small perturbative coupling. One can attempt to solve eg. ([70|) 
in a perturbative expansion in F. That is, set q(t) = go(^) + ^liit) + ■ ■ ■. The solution for 
gi(t) may be found by the Laplace transform: 

gi(t) oc ncos(t) (71) 

clearly exhibiting resonant behavior. This is recognized as a secular term. If eq. (|70D had 
been obtained as an effective equation of motion in a perturbative expansion in A, this would 
be the consistent manner to solve this equation. However, we would be led to conclude that 
perturbation theory breaks down at long times. The correct solution is 

q(t) = e-^*cos[cj(r2)t] ^ cos(t) - ^tcos(t) + 0(T^). (72) 

We see that the first order correction in F is correctly reproduced by perturbation theory, but 
in order to find appreciable damping, we must wait a time ~ C(l/r) at which perturbation 
theory becomes unreliable. 

In order to properly describe dissipation and damping one must resum the perturbative 
expansion. One could in principle keep the first order correction and exponentiate it in 
an ad-hoc manner with the hope that this would be the correct resummation. Although 
ultimately this may he the correct procedure, it is by no means warranted in a field-theoretic 
perturbative expansion, since in field theory, dissipation is related to particle production and 
open multiparticle channels, both very subtle and non-linear mechanisms. Another hint that 
points to a resummation of the perturbative series is provided by the set of equations eqs.(|47|- 
491) . In eq.(^), the classical solution is a periodic function of time of constant amplitude, 
since the classical equation has a conserved energy. As a consequence, the "potential" in 
the equation for the modes (eq.(^)) is a periodic function of time with constant amplitude. 
Thus although the fluctuations react back on the coarse grained field, only the classically 
conserved part of the motion of the coarse grained field enters in the evolution equations of 
the mode functions. This is a result of being consistent with the loop expansion, but clearly 
this approximation is not energy conserving. 

As we will point out in the next section, in an energy conserving scheme the fluctuations 
and amplitudes will grow up to a maximum value and then will always remain bounded at 
all times. 

Thus in summary for this section, we draw the conclusion that perturbation theory is 
not sufficient (without major ad-hoc assumptions) to capture the physics of dissipation and 
damping in real time. A resummation scheme is needed that effectively sums up the whole 
(or partial) perturbative series in a consistent and/or controlled manner, and which provides 
a reliable estimate for the long-time behavior. 
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The next section is devoted to the analytical and numerical study of some of these 
schemes. 
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3 Non-perturbative schemes I: Hartree approximation 



Motivated by the failure of the loop and amplitude expansions, we now proceed to consider 
the equations of motion in some non-perturbative schemes. First we study a single scalar 
model in the time dependent Hartree approximation. After this, we study an 0(N) scalar 
theory in the large N limit. This last case allows us to study the effect of Goldstone bosons 
on the time evolution of the order parameter. 

In a single scalar model described by the Lagrangian density of eq.(|T9|), the Hartree 
approximation is implemented as follows. We again decompose the fields as $^ = (j) + ip^ 
and the Lagrangian density is given by eq . (|2T|) . The Hartree approximation is obtained by 
assuming the factorization (for both ± components) 



(73) 



Translational invariance shows that {{il!^{x,t))'^{x,t)) can only be a function of time, and 
because this is an equal time correlation function, we have that: 



((^+(f,t))^) = ((r(x,t))') = (v^'(t)) 



(74) 



The expectation value will be determined within a self-consistent approximation. After this 
factorization we find: 



where 



(75) 
(76) 



The resulting Hartree equations are obtained by using the tadpole method (?/'^(x, t)) = as 
before. They are given by: 







+ m20 + ^03^^0(V''(t)) 







(27r)3 





(77) 

(78) 
(79) 
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The initial conditions for the mode functions are 

f/+(0) = l ; f/+(0) = -z^fc(O) (80) 

It is clear that the Hartree approximation makes the Lagrangian density quadratic at the 
expense of a self-consistent condition. In the time independent case, this approximation 
sums up all the "daisy" (or "cactus" ) diagrams and leads to a self-consistent gap equation. 

At this stage, we must point out that the Hartree approximation is uncontrolled in this 
single scalar theory. This approximation does, however, become exact in the — > oo limit 
of the 0{N) model which we will discuss in the next section. 



3.1 Renormalization 



We now analyze the renormalization aspects within the Hartree approximation. To study 
the renormalization we need to understand the divergences in this integral 

2 



1) 



(27r)3 2cufc(0) 

The divergences will be determined from the large-k behavior of the mode functions which 
obey the differential equations obtained from (|7^) with the initial conditions (pO|). By a 
WKB-type analysis (see 1 26, 27] for a detailed description), in the k ^ oo limit, we find 

2 

1 1 r _ X _ x _ 1 1 

(82) 



1 
k 



1 



2cufe(0) 

Inserting these results in the integral, it is straightforward to find the divergent terms and 
we find 



d^k 



Kit) 



In 



+ finite 



(83) 



(27r)3 2cjfc(0) 

where A is an upper momentum cutoff and k, a renormalization scale. 

Now, we are in position to specify the renormalization prescription within the Hartree ap- 
proximation. In this approximation, there are no interactions, since the Lagrangian density is 
quadratic. The nonlinearities are encoded in the self-consistency conditions. Because of this, 
there are no counterterms with which to cancel the divergence and the differential equation 
for the mode functions must be finite. Therefore, it leads to the following renormalization 
prescription: 



A 



m 



A, 



B 



A 



m 



A 



R 



R 



R 



(84) 
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where the subscripts B, R refer to the bare and renormahzed quantities respectively and 
is read from eg . (|83D : 

A 

we obtain 



B 



-^A ^ In 



Using this renormahzation prescription eq. 



R 



+ finite 



5) 



ml + 



167r2 



m 



R 



1 



167r2 



In - 



A 



A 



R 



B 



^6) 
^7) 



and 



R 



167r2 



In ^ 







(27r)3 


2cufc(0) 



rA^- 



In 



It is clear that there is no wavefunction renormahzation. This is a consequence of the 
approximation invoked. There is, in fact, no wavefunction renormahzation in either one- 
loop or Hartree approximation for a scalar field theory in three spatial dimensions. 
With an eye towards the numerical analysis, it is more convenient to write 



mt))R = {{rit))R-mo))R) + {^'m 



R 



(89) 



and perform a subtraction at time t = absorbing {iIj^{0))^ into a further finite renormal- 
ization of the mass term (m| + (V'^(0))ij = M^). 

The renormahzed equations that we will solve finally become 



R 



Ut{t) = 





(90) 
(91) 



and 



{mt))R-{^\^)) 
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X 
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+ 



A; 



167r2 



In (^) (0^(t) - 0^(0)) (92) 
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with the initial conditons for Uj^ (t) : 

f/+(0) = 1 ; f7+(0) = -tMO) ; MO) = + Ml + ^02(0) (93) 

It is worth noticing that there is a weak cutoff dependence on the renormahzed equations 
of motion of the order parameter and the mode functions. This is a consequence of the well 
known "triviality" problem of the scalar quartic interaction in four space-time dimensions. 
This has the consequence that for a fixed renormalized coupling the cutoff must be kept 
fixed and finite. The presence of the Landau pole prevents taking the limit of the ultraviolet 
cutoff to infinity while keeping the renormalized coupling fixed. 

This theory is sensible only as a low-energy cutoff effective theory. We then must be 
careful that for a fixed value of Xr, the cutoff must be such that the theory never crosses the 
Landau pole. Thus from a numerical perspective there will always be a cutoff sensitivity in 
the theory. However, for small coupling we expect the cutoff dependence to be rather weak 
(this will be confirmed numerically) provided the cutoff is far away from the Landau pole. 



3.2 Particle Production 

Before we engage ourselves in a numerical integration of the above equations of motion 
we want to address the issue of particle production since it is of great importance for the 
understanding of dissipative processes. 

In what follows, we consider particle production due to the time varying effective mass 
Ai^{t) in eq.(|76D of the quantum field ip for the single scalar model. 



The Lagrangian density for the fluctuations in the Hartree approximation is given by 
eq. (|75[) . Demanding that the order parameter satisfies its equation of motion implies that 
the linear term in 7] in eq. (|75D vanishes. The resulting Lagrangian density is 



m = \ id.^ f - \M\tW + -M) (94) 



To study the issue of the particle production, it proves convenient to pass to Hamiltonian 
density: 

^ = n^^ - £ = ^nj + i (VV')' + \M'{t)i? - (95) 

Here H^ is the canonical momentum conjugate to il). 

In a time dependent background the concept of particle is ambiguous and it must be 
defined with respect to some particular state. Let us consider the Heisenberg fields at t = 
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written as 



n^(x,o) 



E 



(afc(0)+al,(0)) e 



ik-x 



E 



afc(0)-alfc(0))e 



Ak-x 



(96) 



with co'fc(O) as in eq.(^) and VL the spatial volume. The Hamiltonian is diagonalized at t = 
by these creation and destruction operators: 

1- 



/7(0) = E^fe(0) 

k 



4(0)afc(0) + 



(97) 



Thus we define the Hartree-Fock states at t = as the vacuum annihilated by afc(O) 
together with the tower of excitations obtained by applying polynomials in aj.(0) to this 
vacuum state. The Hartree-Fock vacuum state at t = is chosen as the reference state. As 
time passes, particles (as defined with respect to this state) will be produced as a result of 
parametric amplification [p^, |3^ . We should mention that our definition differs from that of 
other authors ^ in that we chose the state at time t = rather than using the adiabatic 
modes (that diagonalize the instantaneous Hamiltonian). 

We define the number density of particles as a function of time as 



where by definition 



d^k Tr 4(0)afe(0)p(t) 



(2vr)= 



Trp(O) 



d^k Tr 4(t)afe(t)p(0) 



(2vr) 



Trp(O) 



4(t) = ?7~^(t,0)4(0)f/(t,0) ; akit) = [/"^(t, 0)afc(0)f/(t, 0) 



(98) 



(99) 



are the time-evolved operators in the Heisenberg picture. In terms of these time evolved 
operators, we may write: 



ipix, t) 



H^(f, t) 
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E 



1 



n T ^/2a;,(0) 



ak{t) + aLfc(t) ) e 



ik-x 



E 



/c^fe(O) 



Ak-x 



(100) 



On the other hand, we now expand the Heisenberg fields at time t in the following 
orthonormal basis 



ipix, t) 
H^(x, t) 



1 
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2cufe(0) 
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k-x 



k-x 



(101) 
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where the mode functions U^{t) ; Uf^{t) = (f/^(t))* are the Hartree-Fock mode functions 
obeying eqs.(91-93) together with the self consistency condition. 

Thus dk, dl are the annihilation and creation operators of Hartree-Fock states, and the 
Heisenberg field ilj{x,t) is a solution of the Heisenberg equations of motion in the Hartree 
approximation. Therefore the al and dk do not depend on time and are identical to a[,(0) 
and ak{0) respectively ( we can check this by evaluating the expansion in eq.( |101| ) at t = 
together with the initial conditions on U'^{t) in eq.(^)). Using the Wronskian properties 
of the function U^{t), we see that the d\ and d^ satisfy the usual canonical commutation 
relations. The reason for the choice of the vacuum state at t = now becomes clear; this is 
the initial time at which the boundary conditions on the modes are determined. The mode 
functions if) ; if) are then identified with positive and negative frequency modes at 
the initial time. 

By comparing the expansion in eq.(0) evaluated at time t with that in eq.( |101| ), we find 
that the creation and annihilation operator at time t can be related to those at time t = 
via a Bogoliubov transformatiom: 

akit) = T+,k{t)ak{0) + ^_,fc(t)a^_fc(0) (102) 

The J^±{t) can be read off in terms of the mode functions U^it) 
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At any time t the expectation value of the number operator for the quanta of tp in each 
k-mode is given by 



A4(t) 

After some algebra, we find 

A4(t) = (2|^+,,(t)|' 



Tr 



Trp{0) 

A4(0) + (|^+,,.(t)|' 



(104) 



(105) 

This result exhibits the contributions from "spontaneous" ( proportional to the initial particle 
occupations) and "induced" ( independent of it) particle production. Since we are analyzing 
the zero temperature case with A/fc(0) = only the induced contribution results. 

Before passing on to the numerical analysis it is important to point out that the Hartree 
approximation is energy conserving. The bare energy density is 



1 • 



S = -<p\t) + V{<P{t)) + 



TrHpiO) 
Trp{0) 



(106) 
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with V{(f)) the classical potential and H is given by eq.(]95|) with 11^ , ip expanded in terms 
of the Hartree mode functions and creation and annihilation operator as in eqs. (|101| , p.01| ). 
Using the equations of motion for (j){t) and the mode functions, and after some lengthy but 
straightforward algebra one finds S = 0. 



3.3 Numerical Analysis 
3.3.1 Unbroken Symmetry Case 

In order to perform a numerical analysis it is necessary to introduce dimensionless quantities 
and it becomes convenient to chose the renormalization point k = Mr. Thus we define 



riit) = m 
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and finally, the equations of motion become 



r] + 1] + 



r]"" + griE{T) = 



+ q' + l + r]\T)+g^{T) 



f/+(0) = 1 ; ^U^{0) = -z^g2 + 1 + ^2(0) 



[/+(r) = 



dr 
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In terms of the dimensionless quantities we obtain the number of particles within a 
correlation volume N{t) = J\f{t)/Ml 



N{t) 



1 pA/Mr 
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111) 



Figures (5.a,b,c) show ^^(r), S(r) and N{t) in the Hartree approximation, for g = 0.1, 
77(0) = 1.0 and A/Mr = 100; we did not detect an appreciable cutoff dependence by varying 
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the cutoff between 50 and 200. Clearly there is no appreciable damping in ^^(t). In fact 
it can be seen that the period of the oscillation is very close to 2tt, which is the period of 
the classical solution of the linear theory. This is understood because the coefficient of the 
cubic term is very small and (7S(r) is negligible. Particle production is also negligible. This 
situation should be contrasted with that shown in figures (6.a-c) and (7.a-c) in which there 
is dissipation and damping in the evolution of //(r) for //(O) = 4, 5 respectively and the 
same values for g and the cutoff. There are several noteworthy features that can be deduced 
from these figures. First the fluctuations become very large, such that (7S(r) becomes 0{1). 
Second, figures (6. a) and (7.a) clearly show that initially, channels are open and energy is 
transferred away from the g = mode of the field. Eventually however, these channels 
shut off, and the dynamics of the expectation value settles into an oscillatory motion. The 
time scale for the shutting off of the dissipative behavior decreases as 7]{0) increases; it is 
about 25 for ri{0) = 4, and about 18 for ri{0) = 5. This time scale is correlated with the 
time scale in which particles are produced by parametric amplification and the quantum 
fluctuations begin to plateau (figures (6.b,c) and (7.b,c)). Clearly the dissipative mechanism 
which damps the motion of the expectation value is particle production. Furthermore, the 
long time dynamics for the expectation value does not correspond to exponential damping. In 
fact, we did not find any appreciable damping for r > 70 in these cases. It is illuminating to 
compare this situation with that of a smaller coupling depicted in figures (8.a-c) for r]{0) = 5 
g = 0.05. Clearly the time scale for damping is much larger and there is still appreciable 
damping at r ^ 100. Notice also that 5'S(r) ^ 2 and that particles are being produced even 
at long times and this again correlates with the evidence that the field expectation value 
shows damped motion at long times, clearly showing that the numerical analysis has not 
reached the asymptotic regime. 

These numerical results provide some interesting observations. First we see that ^^^(t) + 
(7S(r) has asymptotic oscillatory behavior in r. In general this would imply unstable bands 
of wavevectors for the solutions of equation ( |108| ). However, if these unstable bands are 



present, the mode functions would grow and and S(r) would grow as a result since it is an 
integral over all wavevectors. The fact that this does not happen implies that there are no 
unstable bands. 

We have been able to find analytically an exact oscillatory solution of the equation (|108| ) 
which is given by 

l + r]-'{r)+gJ:{T) = -eg - 2 P(r + c^') (112) 



r]{T) = Asn(v/i7^r) (113) 
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where V is the Weierstrass function which is double periodic with periods uo ; uj'^^- The 
parameters in this solution are functions of the coupling. We found that, asymptotically, the 
numerical solution of equations (|108|, |109|) could be fit quite precisely to the above analytic 



solution for precise values of the parameters. These values encode the information about the 
initial conditions and the coupling constant. For this potential given by (|112|) we can find 



the exact solutions for the mode functions in terms of Weierstrass cr functions[^]. In this 
case there are no unstable modes for real wavector q. This is an important result because it 
provides an analytic understanding for the lack of forbidden bands. 

Although this exact result does not illuminate the physical reason why there are no 
unstable bands, we conjecture that this is a consequence of the conservation of energy in the 
Hartree approximation. 

Numerically we find no evidence for exponential relaxation asymptotically, of course, we 
cannot rule out the possibility of power law asymptotic decay for t] with an exceedingly small 
power. 

The fundamental question to be raised at this point is: what is the origin of the damp- 
ing in the evolution of the field expectation value? Clearly this is a coUisionless process as 
collisions are not taken into account in the Hartree approximation (although the one-loop 
diagram that enters in the two particle collision amplitude with the two-particle cut is con- 
tained in the Hartree approximation and is responsible for thresholds to particle production). 
The physical mechanism is very similar to that of Landau damping in the coUisionless Vlasov 



equation for plasmas | TO] and also found in the study of strong electric fields in reference [|n|]. 
In the case under consideration, energy is transfered from the expectation value to the quan- 
tum fluctuations which back-react on the evolution of the field expectation value but out 
of phase. This phase difference between the oscillations of i]^{t) and those of S(r) can be 
clearly seen to be vr in figures (6.a,b), (7.a,b) since the maxima of ^^^(r) occur at the same 
times as the minima of S(r) and viceversa. 

This is an important point learned from our analysis and that is not a priori taken into 
account in the usual arguments for dissipation via collisions. The process of thermaliza- 
tion, however, will necessarily involve collisions and cannot be studied within the schemes 
addressed in this paper. 

3.3.2 Broken Symmetry 

The broken symmetry case is obtained by writing M|. = — < and using the scale 
instead of Mr to define the dimensionless quantities as in eq. ( |107| ) and the renormalization 
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scale. The equations of motion in this case become 
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1 - Mn ( ^ 
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l + ri\r)+gnr) 
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U^{0) = 1 ; ^C/+(0) = -z 



g2 + l + ry2(o) 



(114) 
(115) 



with S(r) given in eg. ( |110| ) but with Mr replaced by nn- The broken symmetry case is more 
subtle because of the possibility of unstable modes for initial conditions in which 77(0) << 1. 
We have kept the boundary conditions of eq.( |115D the same as in eq.( |109| ). This corresponds 
to preparing an initial state as a Gaussian state centered at t]{0) with real and positive 
covariance (width) and letting it evolve for t > in the broken symmetry potential [p6l, ^|. 
The number of particles produced within a correlation volume (now /i/j^) is given by eg. ( |1 1 1| ) 
with Mr replaced by fiR. 

Figures (9.a-c) depict the dynamics for a broken symmetry case in which ri{0) = 10~^ 
i.e. very close to the top of the potential hill. Notice that as the field expectation value rolls 
down the hill the unstable modes make the fluctuation grow dramatically until about r ^ 50 
at which point 5'S(r) ^ 1. At this point, the unstable growth of fluctuations shuts off and 
the field begins damped oscillatory motion around a mean value of about 1.2. This point 
is a minimum of the effective action. The damping of these oscillations is very similar to 
the damping around the origin in the unbroken case. Most of the particle production and 
the largest quantum fluctuations occur when the field expectation value is rolling down the 
region for which there are unstable frequencies for the mode functions (see eq. (|1 14|) ) . This 
behavior is similar to that found previously by some of these author spBj. 



4 Non-perturbative schemes II: Large N limit in the 
0{N) Model 

Although the Hartree approximation offers a non-perturbative resummation of select terms, 
it is not a consistent approximation because there is no a priori small parameter that defines 
the approximation. Furthermore we want to study the effects of dissipation by Goldstone 
bosons in a non perturbative but controlled expansion. 

In this section, we consider the 0{N) model in the large N limit. The large N limit has 
been used in studies of non-equilibrium dynamics p9|, p] , and it provides a very pow- 
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erful tool for studying non-equilibrium dynamics non-perturbatively in a consistent manner. 
The Lagrangian density is the following: 

V{a,^) = imV- 0+^(0-0)' (116) 

for A fixed in the large N limit. Here is an 0{N) vector, = (cr, vr) and vf represents the 
— 1 pions. In what follows, we will consider two different cases of the potential V^((7, tt), 
with (m^ < 0) or without (m^ > 0) symmetry breaking. 
Let us define the fluctuation field xi^yt) 

a{x,t) = (Toit) + x{x,t) (117) 

Expanding the Lagrangian density in eg. ( P-16|) in terms of fluctuations x(2^)^),we obtain 



{(x+-X-),(vf^-vr-)} (118) 



+ 2 



The tadpole condition (x^(x, t)) = will lead to the equations of motion as previously 
discussed. We now introduce a Hartree factorization. In the presence of a non-zero expecta- 
tion value the Hartree factorization is subtle in the case of continuous symmetries. A naive 
Hartree factorization violates the Ward identities related to Goldstone's theorem. This short- 
coming is overcome if the Hartree factorization is implemented in leading order in the large 
N expansionp4||. We will make a series of assumptions that seem to be reasonable, but that 



can only be justified a posteriori when we recognize that with these assumptions, we obtain 
the equations of motion that fulfill the Ward identities. These assumtions are: i) no cross 
correlations between the pions and sigma field and ii) that the two point correlation functions 
of the pion field are diagonal in the 0{N — 1) space given by the remaining unbroken symme- 
try group. Based upon these assumptions we are led to the following Hartree factorization 
of the non-linear terms in the Lagrangian density (again for both ± components): 

^ Q{x^)x^ + constant (119) 
^ 3(x')x (120) 
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(vf )7f + constant 



(121) 

(122) 
(123) 



where by "constant" we mean the operator independent expectation values of the composite 
operators which will not enter into the time evolution equation of the order parameter. 

In this approximation, the resulting Lagrangian density is quadratic, with a linear term 
in X : 

£[ao + x+,^1 - >CK + x-,^1 = {^(9,x+)' + ^(9^^+)'-X^n^) 

- {Mmx^f - \Mlm^f] - {(X^ - X-) , - vf-)} 

(124) 

where V' is the derivative of the Hartree potential defined below. To obtain a large N limit, 
we define 



ao{t) = (l){t)^/N 



(125) 
(126) 



with 

^ 0{1) ; (x') « 0{1) ; « 0(1). (127) 

We will approximate further by neglecting the 0{^) terms in the formal large N hmit. We 
now obtain 
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Using the tadpole method, we obtain the following set of equations: 



; mt)) 

U^{t)^0;u;l,{t)^P + Ml{t) 



d'k U+{t) 



(27r)3 2a;|,(0) 



(128) 

(129) 
(130) 



(131) 
(132) 
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It is clear from the above equations that the Ward identities of Goldstone's theorem are 
indeed fulfilled. Because V'{(p{t),t) = -/N(t){t)M%{t), whenever V'{(t){t),t) vanishes for 7^ 
then A^^ = and the "pions" are the Goldstone bosons. The initial conditions for the mode 
functions U^{t) are 

f/+(0) = l ; f7+(0) = -z^.^fc(O) (133) 

Since in this approximation, the dynamics for the vf and x fields decouples, and the dynamics 
of X does not infiuence that of or the mode functions and {ip"^) we will only concentrate 
on the solution for the vr fields. 



4.1 Renormalization 

The renormalization procedure is exactly the same as that for the Hartree case in the previous 
section (see eqs. ( ^^ - |88|) . We carry out the same renormalization prescription and subtraction 
at t = as in the last section. Thus we find the following equations of motion 



f/+(t)=0 (134) 



with the initial conditions given by eq. ( |133 ) and with the subtracted expectation value given 
byeq.(|92|). 

In contrast with the Hartree equations in the previous section, the cutoff dependence in 
the term proportional to 0^ in eq. (|134|) has disappeared. This is a consequence of the large 
N limit and the Ward identities, which are now obvious at the level of the renormalized 
equations of motion. There is still a very weak cutoff dependence in eq. (|9^) because of the 
triviality issue which is not relieved in the large N limit, but again, this theory only makes 
sense as a low-energy cutoff theory. 



4.2 Numerical Analysis 
4.2.1 Unbroken symmetry 



To solve the evolution equations in eqs. (92, 134) numerically, we now introduce dimensionless 
quantities as in eq. (|107|) obtaining the following dimensionless equations: 

-r^r] + ri + rf + gV^iT) = 
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For particle production in the 0{N) model, the final expression of the expectation value 
of the number operator for each pion field in terms of dimensionless quantities is the same as 
in eq.( |lll| ), but the mode function U^{t) obeys the differential equation in eq. (|108| ) together 
with the self consistent condition. 

Figures (lO.a-c), (ll.a-c) and (12.a-c) show ri{T) ; S(r) ; N{t) for r/(0) = 1; 2 ^ = 0.1; 0.3 
and A/Mji = 100 (although again we did not find cutoff sensitivity). The dynamics is very 
similar to that of the single scalar field in the Hartree approximation, which is not surprising, 
since the equations are very similar (save for the coefficient of the cubic term in the equation 
for the field expectation value). Thus the analysis presented previously for the Hartree 
approximation remains valid in this case. 



4.2.2 Broken symmetry 

The broken symmetry case corresponds to choosing M|, = —fi% < 0. As in the case of the 
Hartree approximation, we now choose fi^ as the scale to define dimensionless variables and 
renormalization scale. The equations of motion for the field expectation value and the mode 
functions now become 



—r]-ri + 7]^ + gT]T.{T) 







^ + q^-l + r]''iT)+gi:ir) 



(r)=0 



(138) 



with S(r) given by eq. ( |137| ). As in the Hartree case, there is a subtlety with the boundary 
conditions for the mode functions because the presence of the instabilities at r = for the 
band of wavevectors < < 1, for ^^^(O) < 1. Following the discussion in the Hartree case 
(broken symmetry) we chose the initial conditions for the mode functions as 



^^(0) = 1 ; |:t/+(0) = -^v'g2 + i + r/2(o) 



(139) 
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These boundary conditions correspond to preparing a gaussian state centered at ri{0) at 
r = and letting this initial state evolve in time in the "broken symmetry potential" . 

Figures (13.a-c) show the dynamics of ?7(t), S(r) and N{t) for ri{0) = 0.5 g = 0.1 and 
cutoff A/fif> = 100. Strong damping behavior is evident, and the time scale of damping is 
correlated with the time scale for growth of S(r) and N{t). The asymptotic value of i]{t) 
and S(t), //(oo) and S(oo) respectively satisfy 

- 1 + r/2(oo) + ^S(cx)) = (140) 

as we confirmed numerically. Thus the mode functions are "massless" describing Goldstone 
bosons. Notice that this value also corresponds to V'{(f)) = in eq. (|128|) . An equilibrium 
self-consistent solution of the equations of motion for the field expectation value and the 
fluctuations is reached for r = oo. 

Figure (13. d) shows the number of particles produced per correlation volume as a function 
of (dimensionless) wavevector g at r = 200. We see that it is strongly peaked at g = 
clearly showing that the particle production mechanism is most efficient for long-wavelength 
Goldstone bosons. Figures (13.e-g) show snapshots of the number of particles as a function 
of wavevector for r = 13, 25, 50 respectively, notice the scale. Clearly at longer times, the 
contributions from q becomes smaller. Figures (14.a-c) show the evolution of i]{t), S(r) 
and N{t) for ri{0) = 10~^, g = 10"^. The initial value of rj is very close to the "top" of the 
potential. Due to the small coupling and the small initial value of ?7(0), the unstable modes 
(those for which < 1) grow for a long time making the fluctuations very large. However 
the fluctuation term S(t) is multiplied by a very small coupling and it has to grow for a long 
time to overcome the instabilities. During this time the field expectation value rolls down 
the potential hill, following a trajectory very close to the classical one. The classical turning 
point of the trajectory beginning very near the top of the potential hill, is close to rjtp = \p2. 
Notice that //(r) exhibits a turning point (maximum) at r/ «i 0.45. Thus the turning point of 
the effective evolution equations is much closer to the origin. This phenomenon shows that 
the effective (non-local) potential is shallower than the classical potential, with the minimum 
moving closer to the origin as a function of time. 

If the energy for the field expectation value was absolutely conserved, the expectation 
value of the scalar field would bounce back to the initial point and oscillate between the 
two classical turning points. However, because the fiuctuations are growing and energy is 
transferred to them from the g = mode, r/ is slowed down as it bounces back, and tends 
to settle at a value very close to the origin (asymptotically about 0.015). Figure (14. b) 
shows that the fluctuations grow initially and stabilize at a value for which gY,{po) ~ 1. 
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The period of explosive growth of the fluctuations is correlated with the strong oscillations 
at the maximum of r){T). This is the time when the fluctuations begin to effectively absorb 
the energy transferred by the field expectation value and when the damping mechanism 
begins to work. Again the asymptotic solution is such that — 1 + 77^(00) + gT,{oo) — 0, and 
the particles produced are indeed Goldstone bosons but the value of the scalar field in the 
broken symmetry minimum is very small (classically it would be rjmin — 1) yet dynamically, 
the field settles at a value 77(00) 0.015!! for g — 10~^). Figure (14.c) shows copious 
particle production, and the asymptotic final state is a highly excited state with a large 
number (C(10^)) of Goldstone bosons per correlation volume. 

The conclusion that we reach from the numerical analysis is that Goldstone bosons are 
extremely effective for dissipation and damping. Most of the initial potential energy of the 
field is converted into particles (Goldstone bosons) and the field expectation value comes to 
rest at long times at a position very close to the origin. 

Notice that the difference with the situation depicted in figures (13.a-c) is in the initial 
conditions and the strength of the coupling. In the case of stronger coupling, the fiuctuations 
grow only for a short time because (i'S(r) becomes 0{1) in short time, dissipation begins 
to act rather rapidly and the expectation value rolls down only for a short span and comes 
to rest at a minimum of the effective action, having transferred all of its potential energy 
difference to produce Goldstone bosons. 

Figures (15.a-c) show a very dramatic picture. In this case 7/(0) = 10^'^, g = 10^^^. 
Now the the field begins very close to the top of the potential hill. This initial condition 
corresponds to a "slow-roll" scenario. The fiuctuations must grow for a long time before 
(7S(r) becomes C(l), during which the field expectation value evolves classically reaching 
the classical turning point at 77^^ = and then bouncing back. But by the time it gets 
near the origin again, the fiuctuations have grown dramatically absorbing most of the energy 
of the field expectation value and completely damping its motion. In this case almost all the 
initial potential energy has been converted into particles. This is a remarkable result. The 
conclusion of this analysis is that the strong dissipation by Goldstone bosons dramatically 
changes the dynamics of the phase transition. For slow-roll initial conditions the scalar 
field relaxes to a final value which is very close to the origin. This is the minimum of 
the effective action, rather than the minimum of the tree-level effective potential. Thus 
dissipative effects by Goldstone bosons introduce a very strong dynamical correction of the 
effective action leading to a very shallow effective potential (the effective action for constant 
field configuration). The condition for this situation to happen is that the period of the 
classical trajectory is of the same order of magnitude as the time scale of growth for the 
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fluctuations. 

The weak coupling estimate for this dynamical non-equilibrium time scale is ~ ln(l/g')/2, 
which is obtained by requiring that 5'S(t) ~ 1. For weak coupling the mode functions grow 
as U^{t) e"^ for < 1, and E(t) e^"^. The numerical analysis confirms this time scale 
for weak coupling. For weakly coupled theories, this non-equilibrium time scale is much 
larger than the static correlation length (in units in which c—1) and the only relevant time 
scale for the dynamics. 

Clearly the outcome of the non-equilibrium evolution will depend on the initial conditions 
of the field expectation value. 

Our results pose a fundamental question: how is it possible to reconcile damping and 
dissipative behavior, as found in this work with time reversal invariance? 

In fact we see no contradiction for the following reason: the dynamics is completely 
determined by the set of equations of motion for the field expectation value and the mode 
functions for the fiuctuations described above. These equations are solved by providing non 
equilibrium initial conditions on the field expectation value, its derivative and the mode 
functions and their derivatives at the initial time t = 0. The problem is then evolved in time 
by solving the coupled second order differential equations. We emphasize the fact that the 
equations are second order in time because these are time reversal invariant. Now consider 
evolving this set of equations up to a positive time to, at which we stop the integration and 
find the value of the field expectation value, its derivative, the value of all the mode functions 
and their derivatives. Because this is a system of differential equations which is second order 
in time, we can take these values at to as initial conditions at this time and evolve backwards 
in time reaching the initial values at t — 0. Notice that doing this involves beginning at a 
time ^0 in excited state with (generally) a large number of particles. The conditions at 
this particular time are such that the energy stored in this excited state is focused in the 
back reaction to the field expectation value that acquires this energy and whose amplitude 
will begin to grow. 

It is at this point where one recognizes the fundamental necessity of the "in- in" formalism 
in which the equations of motion are real and causal. 
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5 Conclusions and Discussion 



We have focused on understanding the first stage of a reheating process, that of dissipation 
in the dynamics of the field expectation value of a scalar field via particle production. Start- 
ing from a scalar field theory with no apparent dissipativc terms in its dynamics, we have 
shown how the evolution of expectation value of the scalar field is affected by the quantum 
fluctuations and particle production resulting in dissipative dynamics. 

Wc started our analysis with a perturbative calculation, both in the amplitude of the 
expectation value of the field and to one-loop order. Analytically and numerically we find 
that dissipative processes cannot be described pcrturbatively. A systematic solution to 
the equation of motion reveals the presence of resonances and secular terms resulting in 
the growth of the corrections to the classical evolution. The perturbative study of the 
0(2) linear sigma model reveals infrared divergences arising from the contribution to the 
dissipative kernels from the Goldstone modes which require non-perturbative resummation. 

A Langevin equation was constructed in an amplitude expansion, it exhibits a generalized 
fluctuation-dissipation theorem with non-Markovian (memory) kernels and colored noise 
correlation functions thus offering a more complex picture of dissipation. 

Motivated by the failure of the perturbative approach, we studied the non-equilibrium 
dynamics in a Hartree approximation both in the symmetric as well as in the broken sym- 
metry case. This approximation clearly exhibits the contribution of open channels and the 
dissipation associated with particle production, which in this approximation is a result of 
parametric ampliflcation of quantum fluctuations. 

In the case of unbroken symmetry we flnd that asymptotically the expectation value of 
the scalar fleld oscillates around the trivial vacuum with an amplitude that depends on the 
coupling and initial conditions. An extensive numerical study of the renormalized equations 
of motion was performed that shows explicitly the dynamics of particle production during 
these oscillations. 

Although the Hartree approximation offers a self-consistent non-perturbative resumma- 
tion, it is not controlled. Thus we were led to study the large N limit in an 0(N) model, 
which also allows us to study in a non-perturbative manner the dissipative dynamics of 
Goldstone bosons. In the case of unbroken symmetry the results are very similar to those 
obtained in the Hartree approximation. 

The broken symmetry case provides several new and remarkable results. 

An extensive numerical study of the equations explicitly shows copious particle produc- 
tion while the expectation value of the field relaxes with strongly damped oscillations towards 
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a minimum of the effective action. 

It is intuitively obvious that Goldstone bosons are extremely effective for dissipation 
since channels are open for arbitrarily small energy transfer. What is remarkable is that 
for "slow-roll" initial conditions, that is the initial expectation value of the scalar field very 
close to the origin and extremely small coupling constant (0(10^^) or smaller) the final value 
of the expectation value remains very close to the origin, most of the potential energy has 
been absorbed in the production of long-wavelength Goldstone bosons. This is confirmed 
by an exhaustive numerical study, including snapshots at different times of the number of 
particles produced for different wavelengths showing a large peak at small wavevectors at long 
times. The numerical study clearly shows that particle production is extremely effective for 
long-wavelength Goldstone bosons. Another remarkable result is that the asymptotic value 
expectation value of the scalar field depends on the initial conditions. These asymptotic 
values correspond to minima of the effective action. Thus we reach the unexpected conclusion 
that in this approximation the effective action has a continuum manifold of minima which 
can be reached from different initial conditions. 

It is pointed out that the basic mechanism for dissipation in this approximation is that 
of Landau damping through the parametric amplification of quantum fluctuations, and pro- 
duction of particles. These fluctuations react-back in the evolution of the expectation value 
of the scalar field, but out of phase. This is a coUisionless mechanism, similar to that found 
in the coUisionless Boltzmann-Vlasov equation for plasmas. 

Our study also reconciles dissipation in the time evolution for the coarse grained variable 
with time reversal invariance, as the evolution is completely specified by an infinite set of 
ordinary second order differential equations in time with proper boundary conditions. 

Our formalism and techniques are sufficiently powerful to give a great deal of insight 
into the particulars of the dissipation process. In particular, we can see that the damping 
of the field expectation value ends as the particle production ends. This shows that our 
interpretation of the damping as being due to particle production is accurate. 

It is useful to compare what we have done here with other work on this issue. We have 
already mentioned the work of Calzetta, Hu|^ and Paz[^. These authors use the closed time 
path formalism to arrive at the effective equations of motion for the expectation value of the 
field. Then they solve the perturbative equations and find dissipative evolution at short times. 
In particular Paz finds the kernel that we have found for the effective equations of motion 
of the expectation value in the perturbative and Hartree case. However, his perturbative 
solution is not consistent. 

We remedy this situation by studing the non-perturbative Hartree equations, which must 
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necesarily be solved numerically as we do. 

There has been other, previous, work on the reheating problem, most notably by Abbott, 
Farhi and Wise[0, Ringwald||45|| and Morikawa and Sasaki In all of these works, the 
standard effective action is used, so that the expectation value is of the "in-out" type and 
hence the equations are non-causal and contain imaginary parts. In essence, they "find" 
dissipational behavior by adding an imaginary part to the frequency that appears in the 
mode equations. We see from our work (as well as that of Calzetta, Hu and Paz) that this 
is not necessary; dissipation can occur even when the system is evolving unitarily. This 
comment deserves a definition of what we call dissipation here: it is the energy transfer from 
the expectation value of the g = mode of the scalar field to the quantum fluctuations 
(g 7^ 0) resulting in damped evolution for the q = mode. 

To what extent are we truly treating the reheating problem of inflationary models? As 
stated in the introduction, reheating typically entails the decay of the infiaton into lighter 
particles during its oscillations. What we do here is understand how the quantum fluctuations 
and the ensuing particle production influence the dynamics of the evolution of the expectation 
value of the field. Thus technically speaking, this is not the reheating problem. However, 
we are able to understand where dissipation comes from in a field theory, and are able to 
give a quantitative description of the damping process for the expectation value of the field. 
We expect that the physics of dissipation when the scalar field couple to others will be very 
similar to the case studied in this article. 

Furthermore, the techniques we develop here are easily adapted to the case where the 
infiaton couples to fermions, a case that we are currently studying|^6[. In this connection, 
during the writing of this paper, two related pieces of work on the reheating problem have 
appeared[^]. They both look at the effect of particle production from the oscillations of 
the infiaton field due to parametric amplification. What they do not do is to account for 
the back reaction of the produced particles on the evolution of the expectation value of the 
infiaton. As we have learned with our study, this back reaction will eventually shut off the 
particle production, so that these authors may have overestimated the amount or particle 
production. 

We recognize that our non-perturbative treatment neglects the effect of collisions as 
mentioned above. Dissipation appears in a manner similar to Landau damping. In the 
Hartree approximation there does not seem to be a natural way to incorporate scattering 
processes because this is a mean-field theory. However the large N expansion allows a 
consistent treatment of scattering processes for which the first contribution (2-2 particle 
scattering) will appear at 0{1/N). This will be a necessary next step in order to fully 
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understand the coUisional thermalization which is the second stage of the reheating process. 
This will clearly be a fascinating and worthy endeavour that we expect to undertake soon. 
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7 Appendix I 



In this appendix we give an alternative but equivalent method to derive the equations of 
motion based on the direct time evolution of initial density matrices. This method gives rise 
to equations identical in form to those given in section 2 and 3. However, this method might 
have greater applicability for initial density matrices that are not of the thermal type, such 
as might appear in chaotic inflation models, but is restricted to the gaussian approximation. 



For more details the reader is referred tolET 



Our starting point is again the Liouville equation for the density matrix of the system: 

Z^^ = [iJevol,p(t)], (141) 

where Hemi is the Hamiltonian of the system that drives the non- equilibrium evolution of 
the system. 

Next we use the density matrix to extract the order parameter from the field in the 
Schroedinger picture: 

0W = ^/ d'xTilpmix)], (142) 

with Q the spatial volume (later taken to infinity), and $(a;) the field in the Schroedinger 
picture. Using the Liouville equation together with the Hamiltonian: 

H = /rf^x|H! + l(V<l>)^ + n<l>)| (143) 

and the standard equal time commutation relations between a field and its canonically con- 
jugate momentum, we find the equations: 

^ = i/A{n(f,«)) = i/Arr|p(«)n(f)i = i(i) (lu) 

= 4 f,M'-^). (145) 
From these equations we can find the equation of motion for the order parameter 0(t): 

We expand the field operator as $(x) = (pit) +ip{x, t), insert this expansion into equation 
( |146| ) and keep only the quadratic terms in the fluctuation field '?/'(x, t). Doing this yields 
the equation: 
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^ + nm) + / , 0) + . . . = (147) 

Here the primes stand for derivatives with respect to (p. 

To make sense of the above equation, we need to specify how the expectation value is to 
be taken. This entails a specification of the basis in functional space we will use to write the 
density matrix in, and, once given this basis, what the form of the density matrix will be. 

In order to generate an expansion in h, we choose the density matrix to be a Gaussian 
with a covariance of order 1/h, in the basis given by that of the spatial Fourier components 
of the fluctuation field ■ip{x,t). Thus write ip{x,t) as: 

= ^T.Mt)e-''-' (148) 

If we are choosing the state to be Gaussian, we need to make the Hamiltonian quadratic, in 
order that the time evolved state remain Gaussian, thus allowing for a consistent approxi- 
mation scheme. To do this we take the Hamiltonian above and expand it in terms of the 
fluctuation field out to quadratic order. If we then insert the expansion in equation ( |148[ ) 
into the truncated Hamiltonian obtained in this way, we arrive at: 



H. 



quad 



f!V(,#,(i)) + 



2 4^ 

k 



6^ 



2vmm-k+^iit)i^k^- 



ulit) = p + v"{m)- 

For the density matrix in this basis, we make the following ansatz\ 



(149) 



(150) 



^^,{t)^_,{t) + ^^,(t)^_,(t) + 



Bk{t) 



Mi) 



(151) 

(152) 
(153) 



where 0(t) = ($(a;)) and 7Tk{t) is the Fourier transform of (H(x)). This form of the density 
matrix is dictated by the hermiticity condition p'''[$, = p*[$,$,t]; as a result of this. 



46 



Bk{t) is real. The kernel Bk{t) determines the amount of "mixing" in the density matrix, 
since if 1?^ = 0, the density matrix describes a pure state because it is a wave functional 
times its complex conjugate. 

In the Schrodinger picture, the Liouville equation ( |141j ) becomes 

+V'.k{<P{t)) {^Pk - ^k) + ^ujlit) (^,^^_,^ - i^k^P^k) } p[$, t] (154) 

Since the modes do not mix in this approximation to the Hamiltonian, the equations for 
the kernels in the density matrix are obtained by comparing the powers of ip on both sides 
of the above equation. We obtain the following equations for the coefficients: 



.A 



Al-Bl-ujlit) 
BkiAk - Al) 



TC 



(155) 

(156) 
(157) 
(158) 
(159) 



The last two equations are identified with the classical equations of motion for the order 
parameter ( [147| ). The equation for Bk{t) refiects the fact that a pure state Bk = remains 
pure under time evolution. 

At this point we have not specified any initial conditions; if we were to specify that the 
initial state is thermal, we would be led to the following initial conditions for the kernels 
Ak{t), Bk{t), as well as for the order parameter and its momentum: 



Ak{0) 
BkiO) 

A4(0) 

m 



sinh [(3ohwk{0)] 



tanh 



TTh 

(j)o ; 7r(0) = TTo 



Po^WkiO) 



(160) 
(161) 

(162) 
(163) 



Here (3o is the inverse temperature at the initial time (which we have taken to infinity in this 
work), while the frequencies Wk{0) are given by: 



47 



wliO) = e + V"iMto)) (164) 

Note that up to this point, the formahsm we have derived is quite general; it is only in 
the specification of the initial conditions that the thermal density matrix comes in. This 
is to be contrasted with the closed time path method, which implicitly assumes a thermal 
density matrix for the initial state. 

Defining a complex function 

Ak{t) = AkR{t) + Aki{t) 
AkRit) = AkRit) tanhPoh WkiO) = -Bk{t) smhpoTi Wk{Q) 
Aki{t) = Akiit) 
Afl(O) = Wk{^) 

A/(0) = 0, (165) 
we find that Akif) satisfies the following Ricatti equation: 

lAkit) = Akitf - Wk{tf (166) 

This equation can be recast in a form closer to that for the mode functions of section 2 
by the following set of transformations: 



lAkit) ^ ^ 
Xk 

Xk{t) = I ^ {Uik + OA^k) (167) 



From the thermal initial conditions, and the equation for Ak{t), we find that the mode 
functions Uak satisfy the following equations: 

{^ + ^fc(i)Kfc(t) = 

Wifc(O) = 1, ; Wifc(O) = 

W2fc(0) = 0, ; U2k{Q)=Wk{Q) (168) 

Finally, using these mode functions to calculate the expectation value of il){x,tY that 
enters into the equation of motion for the field expectation value, we find that equation 
( |146| ) becomes: 
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cj>{t) + V'{cj>{t)) + 



2 2 J (27r)3 Wfc(O) 



coth 



2 



= 



(169) 



In the limit that the initial temperature goes to zero, these equations are the same as 
the one-loop equations derived in section 2. To obtain the evolutions equations in the 
Hartree approximation or in the large N limit in the 0{N) theory one invokes the Hartree 
factorization and proposes a gaussian ansatz as in the one loop case solving for the covariance 
following the same steps as above. We have included this derivation since it may be applicable 
to more general situations, and to give an alternative to the complex path formulation for 
the gaussian approximation in non-equilibrium problems . 



49 



8 Appendix II 



In this appendix we provide some necessary details that may help the reader reproduce our 



appearing in various models under different approximation schemes. By examing the cor- 
responding differential equations (for example, eq. ( |135| )) with the corresponding boundary 
conditions a typical integration routine will face extremely large derivatives when the values 
of the wavevectors q reach the cutoff. We have used improved fourth order Runge-Kutta 
algorithms for the numerical integration, but any algorithm will face the same problems for 
large enough cutoff. We recognize however, that for large q (typically of the order of the 
dimensionless cutoff momentum ^ 100) for which q^ » 1 + if'{T) + g^i^) the solutions to 
these equations are like plane waves with fast varying phases in time. It is precisely this fast 
variation of the phase that leads to potential numerical problems. We introduce a WKB-like 
ansatz to take out the fast varying phases from U^{t) by defining new mode functions which 
are smoothly varying in time for large q. We consider specifically the 0{N) model in the 
large N limit as an example. It is straightforward to apply this analysis to the other cases 
in this paper. 

First of all, we define the new mode functions U^{t) as follows: 



numerical results. These concern the numerical integration for the mode functions U^{t) 



W+(r) = e^/o^'^-'^«(-')f/;(r) ; Cu, = ^q^ + 1 + v^r) (170) 



The initial conditions for U^{t) defined above can be read off from (|115|) for the symmetric 



phase: 




the slowly varying functions obey: 



In can be seen that in the large q limit, (g^ >> 1 + i]'^{t) + g^ir) 



(172) 



(r) ^ ; for T>0 



(173) 



which are smooth functions at all times for large q. 

For the broken symmetry case in 0{N) model, the mode functions U^{t) are solutions 
of the differential equation 



iu,iT)U^ + [-2 + g^{T)]U^ = 



(174) 
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with the same initial conditions as in ( |171| ). 
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Figure Captions 

Fig.l: Diagrams contributing to the equation of motion up to one-loop order, and in the 
amphtude expansion up to 

Fig. 2: First order quantum correction for discrete symmetry case 771 (r) for 771(0) = 
; 7)1(0) = ; rjciiO) = 1 ; r)ci{0) = 0. The cutoff is A/ruR = 100. 

Fig. 3: First order quantum correction for 0(2) model 771 (r) for 771(0) = ; 7)1(0) = 
; r]ci{0) = 0.6 ; 7)^(0) = 0. The cutoff is A/hr = 100. 

Fig. 4: One loop diagrams contributing to the effective action up to O((0='=)^). The 
dashed external legs correspond to the zero mode. 

Fig. 5. a: ?7(r) vs r in the Hartree approximation, unbroken symmetry case, g — 
0.1 ; ri{0) = 1 ; A/Mr = 100. 

Fig. 5.b: E(r) vs r for the same case as in fig. (5. a). 

Fig. 5.c: N{t) vs r for the same case as in fig. (5. a). 

Fig. 6. a: r7(r) vs r in the Hartree approximation, unbroken symmetry case, g = 
0.1 ; 7/(0) = 4 ; A/ Mr = 100. 

Fig. 6.b: S(r) vs r for the same case as in fig. (6. a). 
Fig. 6.c: N{t) vs r for the same case as in fig. (6. a). 

Fig. 7. a: ^^(r) vs r in the Hartree approximation, unbroken symmetry case, g = 
0.1 ; 77(0) = 5 ; A/ Mr = 100. 

Fig. 7.b: S(r) vs r for the same case as in fig. (7. a). 
Fig. 7.c: N{t) vs r for the same case as in fig. (7.a). 

Fig. 8. a: 77(r) vs r in the Hartree approximation, unbroken symmetry case, g — 
0.05 ; 77(0) = 5 ; A/Mr = 100. 

Fig. 8.b: S(r) vs r for the same case as in fig. (8. a). 
Fig. 8.c: N{t) vs t for the same case as in fig. (8. a). 

Fig. 9. a: 77(t) vs r in the Hartree approximation, broken symmetry case, g — 10^^ ; 77(0) = 
10-5 ; A/iiR = 100. 

Fig. 9.b: E(r) vs r for the same case as in fig. (9. a). 
Fig. 9.c: N{t) vs t for the same case as in fig. (9. a). 

Fig. 10. a: 77(t) vs t in the large N approximation in the 0(N) model, unbroken symmetry 
case. ^ = 0.1 ; 77(0) = 1 ; A/Mr = 100. 

Fig. lO.b: E(r) vs r for the same case as in fig. (10. a). 
Fig. lO.c: N{t) vs r for the same case as in fig. (10. a). 

Fig. 11. a: 77(t) vs t in the large N approximation in the 0(N) model, unbroken symmetry 
case. ^ = 0.1 ; 77(0) = 2 ; A/Mr = 100. 
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Fig. ll.b: E(t) vs t for the same case as in fig. (11. a). 
Fig. ll.c: N{t) vs t for the same case as in fig. (11. a). 

Fig. 12. a: 77 (r) vs r in the large N approximation in the 0(N) model, unbroken symmetry 
case. ^ = 0.3 ; 77(0) = 1 ; A/Mr = 100. 

Fig. 12.b: E(r) vs r for the same case as in fig. (12. a). 
Fig. 12.c: N{t) vs r for the same case as in fig. (12. a). 

Fig. 13. a: 77(t) vs r in the large N approximation in the 0(N) model, broken symmetry 
case. ^ = 0.1 ; r/(0) = 0.5 ; A/fiR = 100. 

Fig. 13. b: S(r) vs r for the same case as in fig. (13. a). 
Fig. 13. c: N{r) vs r for the same case as in fig. (13. a). 

Fig. 13. d: Number of particles in (dimensionless) wavevector q, Nq{T) at r = 200. 
Fig. 13. e: Number of particles in (dimensionless) wavevector q, Nq{r) at r = 13. 
Fig. 13. f: Number of particles in (dimensionless) wavevector q, Nq{r) at r = 25. 
Fig. 13. g: Number of particles in (dimensionless) wavevector q, Nq{r) at r = 50. 
Fig. 14. a: 77 (r) vs r in the large N approximation in the 0(N) model, broken symmetry 
case, g = 10"^ ; 77(0) = lO'^ ; A/fiR = 100. 

Fig. 14. b: S(r) vs r for the same case as in fig. (14. a). 
Fig. 14. c: N{t) vs r for the same case as in fig. (14. a). 

Fig. 15. a: 77(r) vs r in the large N approximation in the 0(N) model, broken symmetry 
case, g = IQ-^^ ; 7/(0) = 10"^ ; A/fiR = 100. 

Fig. 15.b: E(r) vs r for the same case as in fig. (15. a). 
Fig. 15.c: N{t) vs r for the same case as in fig. (15. a). 
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